/**
 * Perform cubic spline interpolation to estimate the value at a given point.
 *
 * @param {number} interpolationPoint - The point at which to perform interpolation.
 * @param {number[]} xValues - The array of x values.
 * @param {number[]} yValues - The array of corresponding y values.
 * @returns {number | null} - The interpolated value at the given point or null if outside the data range.
 */
const cubicSplineInterpolation = (interpolationPoint, xValues, yValues) => {
    const numPoints = xValues.length - 1
    const intervals = []
    const slope = []
    const lambda = []
    const mu = []
    const zeta = []
    const coefficientsC = []
    const coefficientsB = []
    const coefficientsD = []

    for (let i = 0; i < numPoints; i++) {
        intervals[i] = xValues[i + 1] - xValues[i]
        slope[i] = (3 / intervals[i]) * (yValues[i + 1] - yValues[i]) - (3 / intervals[i - 1]) * (yValues[i] - yValues[i - 1])
    }

    lambda[0] = 1
    mu[0] = 0
    zeta[0] = 0

    for (let i = 1; i < numPoints; i++) {
        lambda[i] = 2 * (xValues[i + 1] - xValues[i - 1]) - intervals[i - 1] * mu[i - 1]
        mu[i] = intervals[i] / lambda[i]
        zeta[i] = (slope[i] - intervals[i - 1] * zeta[i - 1]) / lambda[i]
    }

    lambda[numPoints] = 1
    zeta[numPoints] = 0
    coefficientsC[numPoints] = 0

    for (let j = numPoints - 1; j >= 0; j--) {
        coefficientsC[j] = zeta[j] - mu[j] * coefficientsC[j + 1]
        coefficientsB[j] =
            (yValues[j + 1] - yValues[j]) / intervals[j] -
            (intervals[j] * (coefficientsC[j + 1] + 2 * coefficientsC[j])) / 3
        coefficientsD[j] = (coefficientsC[j + 1] - coefficientsC[j]) / (3 * intervals[j])
    }

    for (let j = 0; j < numPoints; j++) {
        if (interpolationPoint >= xValues[j] && interpolationPoint <= xValues[j + 1]) {
            const deltaX = interpolationPoint - xValues[j]
            return (
                yValues[j] + coefficientsB[j] * deltaX + coefficientsC[j] * deltaX ** 2 + coefficientsD[j] * deltaX ** 3
            )
        }
    }

    return null
}


export { cubicSplineInterpolation }
